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FIRE  AND  SMOKE  SIMULATOR  (FSSIM) 
VERSION  1  -  THEORY  MANUAL 


1.1  Background 

In  both  peacetime  and  during  war,  fire  represents  a  significant  threat  to  any  ship.  A  fire, 
whether  started  by  a  mechanical  failure,  intentional  or  unintentional  human  actions,  or  damage 
fi'om  a  weapon  hit,  threatens  the  ship  in  a  number  of  ways.  The  crew’s  health  and  ability  to 
operate  the  ship  is  affected  by  direct  exposure  to  the  fire  or  by  the  spread  of  smoke  and  toxic 
gasses  through  the  ship  by  either  natural  or  mechanical  ventilation.  Electrical  systems  can  be 
degraded  by  thermal  exposure,  exposure  to  acid  gasses  in  the  combustion  products,  or  by 
electrical  failure  resulting  firom  soot  deposition,  which  might  include  hampered  cooling  or 
dielectric  breakdown  fi'om  the  electrical  conductivity  of  the  soot.  Mechanical  systems  can  suffer 
thermal  damage.  Lastly,  fire  growth  and  spread  could  potentially  ignite  explosive  materials, 
rocket  motors,  aviation  fuel,  or  other  highly  flammable  substances,  which  could  possibly  result 
in  temperatures  or  overpressures  high  enough  to  affect  the  ship  structurally. 

There  currently  exist  a  number  of  computational  tools  for  examining  the  effects  of  a  fire  that 
can  be  applied  to  a  ship  and  its  crew.  One  could  use  hand  calculations  for  examining  simple 
scenarios  in  single  compartments.  Simple  rules  can  be  used  to  extend  this  approach  to  multiple 
compartments.  Zone  models  are  suitable  for  examining  somewhat  more  complex, 
time-dependent  scenarios  involving  multiple  compartments  and  levels,  but  stability  can  be  a 
problem  for  multi-level  scenarios,  scenarios  with  Heating,  Ventilation,  and  Air  Conditioning 
(HVAC)  systems,  and  for  post-flashover  conditions.  Computational  fluid  dynamics  (CFD) 
models  can  yield  detailed  information  about  temperatures,  heat  fluxes,  and  species 
concentrations;  however,  the  time  penalty  of  this  approach  currently  makes  using  CFD 
unfeasible  for  long  periods  of  real  time  or  large  computational  domains.  There  is  another  class 
of  models  which  have  traditionally  played  supporting  roles  in  fire  modeling.  This  class  is 
constituted  by  a  variety  of  network  models,  which  are  used  for  ventilation  systems  in  buildings 
or  fluid  flow  in  piping  networks.  These  models,  however,  lack  specific  physics  required  for  fire 
modeling. 

The  need  for  fire  modeling  occurs  through  the  “cradle-to-grave”  lifetime  of  a  ship,  which 
encompasses  the  phases  of  concept  design,  detailed  design  and  hazards  evaluation,  doctrine  and 
tactics  development,  ship  operation  and  maintenance,  and  crew  training.  In  the  concept  design 
phase,  fire  modeling  is  needed  to  evaluate  ship  designs  and  design  philosophies  in  order  to 
quickly  arrive  at  an  overall  concept  to  meet  required  performance  goals.  As  this  concept  is 
refined  into  a  detailed  design,  fire  modeling  is  continued  to  evaluate  ship  vulnerability  and 
recoverability  and  to  begin  the  process  of  defining  ship  operation.  Both  of  these  phases  could 
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potentially  require  evaluating  the  fire  behavior  of  a  ship  and  its  system  for  hundreds  of  scenarios 
in  design  cycles  spanning  a  few  months  time.  Doctrine  and  tactics  planning  also  require 
evaluating  many  separate  fire  scenarios  in  combination  with  a  number  of  candidate  doctrines. 
During  the  operational  phase  of  a  ship’s  lifetime  fire  modeling  can  support  crew  training  through 
interactive  training  simulations.  Fire  modeling  can  also  aid  in  the  recovery  from  a  fire  event  by 
providing  a  faster  than  real  time  predictive  capability  to  gauge  the  potential  effectiveness  of  a 
course  of  action. 

For  each  of  the  aforementioned  uses  of  fire  modeling,  speed  and  the  ability  to  model  large 
complexly  interconnected  spaces  with  ventilation,  detection,  and  suppression  with  a  reasonable 
degree  of  accuracy  are  performance  requirements  for  the  model.  Hand  calculations,  while  fast, 
have  limitations  in  applicability  and  large  imcertainties  in  their  results.  CFD  computations  have 
the  potential  to  be  accurate  to  any  desired  degree,  but  are  slow.  Zone  models  balance  well 
computational  time  and  accuracy;  however,  in  general  these  models  do  not  perform  well  with 
ventilation  systems  and  they  lack  the  ability  to  model  control  systems.  To  meet  the 
computational  speed  and  algorithm  requirements,  it  was  decided  to  develop  Fire  and  Smoke 
Simulator  (FSSIM)  as  a  network  fire  model. 

1.2  Network  Model  Overview 

A  network  model  represents  a  one-dimensional  abstraction  of  the  physical  world.  The 
prototype  being  modeled  is  represented  by  a  set  of  nodes  and  node  connections,  where  the  nodes 
contain  one  set  of  physical  variables  and  the  node  connections  represent  transfers  of  those 
variables  between  nodes.  Contrast  this  with  a  zone  model  which  uses  two  sets  of  variables  or  a 
CFD  model  which  uses  a  multitude.  A  network  representation  allows  for  maximum  physical 
extent  of  a  simulation  with  a  minimum  set  of  equations,  one  per  geometric  volmne  of  interest. 
Since  the  number  of  locations  being  solved  for  is  kept  to  a  minimum,  a  network  model  also  has 
the  potential  for  the  faster  execution  speeds. 

There  are  currently  a  number  of  network  models  in  existence  that  are  used  for  safety  related 
simulations.  Some  of  these  include: 

•  TRAC  -  Transient  Reactor  Analysis  Code.  Used  for  simulating  normal  and  accident 
conditions  of  nuclear  power  plant  cooling  systems  (Spore,  et  al.,  2000). 

•  RELAP  -  Reactor  Excursion  and  Leak  Analysis  Program.  Used  for  simulating 
normal  and  accident  conditions  of  nuclear  power  plant  cooling  systems  (USNRC, 
2001). 

•  GOTHIC  -  Generation  of  Thermal  Hydraulic  Information  for  Containments.  Used 
for  simulating  normal  and  accident  conditions  of  nuclear  power  plant  containment 
buildings  and  reactor  buildings  (George,  et  al.,  2000a-b). 

•  MELCOR  -  Used  for  simulating  normal  and  accident  conditions  of  nuclear  power 
plant  containment  buildings  and  reactor  buildings  including  aerosol  dispersion,  core 
melt/core  melt  interactions,  and  hydrogen  combustion  (Gauntt,  et  al.,  2000). 

•  CONTAM  -  Used  for  simulating  the  dispersal  of  airborne  contaminants  through  a 
building  and  its  associated  HVAC  systems  (Dols,  Walton,  and  Denton,  2000). 
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•  COCOSYS  -  containment  COdeSYStems  Used  for  simulating  accident  conditions 
in  light  water  nuclear  power  plant  containment  buildings  and  reactor  buildings 
(Allelein  et  al.,  1999;  Klein-HejSling  et  ah,  2000). 

1.3  FSSIMvl.O  Overview 

FSSIM  is  a  network  fire  model  written  to  simulate  the  spread  of  fire  and  smoke  in  a  naval 
vessel.  However,  there  is  nothing  in  the  model  to  preclude  its  use  for  other  types  of 
compartmented  structures  such  as  a  building.  FSSIM  is  written  in  standard  Fortran  95  (Adams, 
et  al.,  1997)  and  as  such  is  capable  of  being  compiled  on  any  platform  for  which  a  Fortran  95 
compiler  is  available.  In  FSSIM  each  compartment  in  a  structure  is  represented  as  a  single  node 
with  surfaces  (e.g.  bulkheads,  decks,  and  overheads)  and  vent  openings  (e.g.  doors,  hatches,  etc.) 
represented  as  node  cormections.  There  is  no  practical  limit  (2x10^  items  of  any  type)  to  the  size 
of  an  FSSIM  simulation  other  than  available  computational  resources  and  time. 

FSSIMvl.O  has  the  following  capabilities: 

•  ID  flow  model  including  fiiction  losses  and  temperature-dependent  specific  heat. 

•  ID  multiple-layer,  temperature-dependent  heat  transfer. 

•  N-surface,  gray-gas  radiation  heat  transfer,  including  radiation  streaming  through 
openings. 

•  Bidirectional  flow  through  horizontal  (hatches)  and  vertical  (doors)  flow  connections. 

•  Combustion  product  species  tracking. 

•  Oxygen  and  fuel  limited  combustion. 

•  Multiple  user-defined  fires  along  with  fire  spread  via  compartment-to-compartment 
heat  transfer. 

•  HVAC  systems  including  ducts,  dampers,  and  fans  with  forward  and  reverse  flow 
losses  and  multiple  fan  models. 

•  Fire  detection  via  heat,  smoke,  and  fire  detection 

•  Fire  spread  by  compartment-specific  criteria. 

•  Fire  suppression  via  sprinklers,  water  mist,  gaseous  agents,  aerosol  agents,  and  foam. 

•  Fire  spread  prevention  via  boundary  cooling. 

•  Simple  control  systems  to  link  operation  of  equipment  to  sensors  or  times. 

•  Fast,  near-real-time  execution  speed. 

1.4  Documentation  Overview 

FSSIM  documentation  consists  of  three  volumes.  The  FSSIM  Theory  Manual  (this 
document),  the  FSSIM  Users’  Manual  (Floyd,  et  al.,  2004),  and  the  FSSIM  Validation  Manual 
(Floyd,  et  al.,  2004).  The  Theory  Manual  describes  the  equations  solved,  and  the  solution 
algorithm  for  the  heat  and  mass  transfer  along  with  the  equations  and  algorithms  for  FSSIM 
submodels.  The  Users’  manual  discusses  creation  of  FSSIM  input  files  including  guidance  in 
transforming  a  prototype  structure  to  a  network  representation  and  explains  the  available  outputs. 
Lastly,  the  Validation  Manual  documents  FSSIM  performance  as  compared  to  experimental  data 
and  other  fire  models. 
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There  is  a  companion  viewer  for  FSSIM  called  [MSU  GUI  NAME]  that  is  under 
development  by  Mississippi  State  University.  This  viewer  will  display  animated  results  of 
FSSIM  calculations  in  a  3D  representation  of  the  ship  (Haupt,  et  al.,  2004). Hydraulic  Solver 

2.0  HYDRAULIC  SOLVER 

2.1  Overview 

FSSIM  uses  a  network  or  lumped-parameter  approach  to  solve  the  time-dependent  flow  of 
gases  through  a  series  of  vent-connected  compartments  and  associated  ventilation  systems.  In 
this  approach,  compartments  are  considered  one  dimensional,  that  is  the  pressure,  mass,  and 
energy  content  of  a  compartment  are  represented  by  one  set  of  integral  quantities.  Flows  through 
vents  can  be  bi-directional  in  the  sense  of  simultaneous  inflow  and  outflow,  but  vector  flow 
quantities,  such  as  absolute  flow  direction,  are  not  tracked.  Ventilation  systems  are  broken  down 
in  the  simplest  possible  duct-node  representation.  Flow  in  ducts  is  assumed  to  be  imidirectional. 
A  duct  is  any  length  of  ducting  with  bends,  expansions,  contractions,  etc.,  where  only  one  flow 
possibility  exists.  A  node  is  typically  a  duct  terminal,  which  can  be  either  a  compartment  or  a 
ventilation  component  joining  multiple  ducts  such  as  a  manifold  or  a  tee.  Figure  1  shows  a 
simple  compartment  layout  and  its  network  representation. 


/  ►  / 

7  1  / 

■m 

imi 

1^11 

hhi 

i^&i 

One  fen 


Five  heat  transfer  Four  heat  transfer  Five  heat  transfer 
connections  to  the  connections  to  the  connections  to  the 
ambient  ambient  ambient 


Figure  1 :  Network  Model  Concept 
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2.2  Conservation  Equations  in  Derivative  Form 

The  FSSIM  hydraulic  solver  solves  the  ID  conservation  equations  for  mass,  momentum  and 
energy.  Energy  and  mass  are  conserved  explicitly,  whereas  momentum  is  conserved  implicitly. 
Energy  and  mass  conservation  use  a  control-volume  approach,  where  the  control  volume  is  either 
a  single  compartment  or  a  ventilation  system  node.  Momentum  is  implicitly  solved  for  at  vent 
connections  or  in  ducts.  The  set  of  equations  used  is  taken  from  MELCOR  (Gauntt,  et  al.,  2000) 
with  the  removal  of  multiphase  (condensable  vapor,  droplet,  liquid  pools,  and  bubbles)  and 
interphase  (phase-related  drag  and  phase  change)  terms. 


Mass  Equation: 


dM, 


dt 


-  =  E'^ijP?VjFA+Mj 


(1) 


In  this  equation  i  is  the  control  volume  (compartment  or  duct  node)  j  is  a  flow  connection,  d 
is  a  donor  or  upstream  quantity,  M  is  the  total  mass  in  a  voliime,  A  is  flow  area,  F  is  a  function 
that  modifies  the  flow  area  (allows  for  opening/closing),  a  is  a  flow  direction  and  is  +1  for  flow 
out  of  i  and  -1  for  flow  into  i,  v  is  the  velocity,  and  p  is  the  density.  The  mass  source  term,  M; , 
represents  mass  additions  due  to  pyrolysis  and  suppression.  The  change  in  mass  over  time  is 
equal  to  the  net  flow  into  the  volume  plus  the  net  source/sink  of  mass  in  the  volume.  For  a  duct 
node,  there  is  no  mass  storage  or  source  term,  so  a  modified  equation  can  be  given: 

0  =  2]^ijPjVA  (2) 


Energy  Equation: 


dE: 


dt 


-  =  Z<>s'’jF,AjpX+E, 


(3) 


In  this  equation,  variable  definitions  are  as  above  with  the  addition  of  E  for  the  energy 
content  of  the  volume,  and  h  is  a  species  enthalpy.  As  with  mass,  duct  nodes  do  not  store 
energy.  The  energy  source  term,  Ej  ,represents  the  addition  or  removal  of  energy  due  to 
combustion  heat  release,  heat  removal  due  to  suppression,  heat  transfer  to  surfaces,  and  enthalpy 
changes  due  to  mass  changes  resulting  from  pyrolysis,  suppression,  and  compartment  vent  flows. 

Momentum  Equation: 


PjL,  ^  =  (P,  -  Pj  +  (pgAz),  +  APj  -  1k,Pj 


V-  V. 

J  J 


(4) 


Additional  variables  in  this  equation  are  k  for  a  control  volume,  P  for  a  control  volume 
pressure,  AP  for  a  momentum  source  term  such  as  a  fan,  K  for  flow  losses  resulting  from  wall 
friction  and  form  losses  (e.g.  fittings,  flow  obstructions,  area  changes,  and  entrances/exits). 
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(pgAz)  for  buoyancy  head,  and  L  for  the  momentum  length  of  the  flow  connection.  Buoyancy 
head  and  momentum  length  will  be  discussed  later.  Note  that  the  above  momentum  equation 
does  not  contain  a  term  for  kinetic  energy.  This  term  is  impractical  to  define  for  a  lumped 
parameter  model,  and  in  general  its  omission  is  a  minor  source  of  error  (Gauntt,  et  al.,  2000) 

In  general,  this  set  of  equations  must  be  solved  in  an  iterative  manner  due  to  the  stif&iess  of 
the  pressure-velocity  coupling.  There  may  be  flow  connections  whose  velocities  oscillate 
direction  during  this  iterative  process.  To  aid  in  convergence  an  iteration-dependent  donor 
formulation  is  used.  In  the  momentum  equation  the  density  term  is  modified  as  follows: 

pr'=pp;+o-pk  (5) 

The  density  used  in  the  momentum  equation  at  inner  loop  iteration  n+1,  see  2.6.2,  is  a  linear 
combination  of  the  density  used  for  iteration  n  and  the  donor  or  upstream  density.  P  ranges  from 
0  to  1 .  For  the  first  third  of  the  iterations  p=0  and  a  pure  upstream  density  is  used.  For  the  last 
third  of  the  iterations  p  =1  and  the  prior  iteration  density  is  used,  p  varies  linearly  from  0  to  1 
over  the  middle  third  of  the  iterations.  If  flow  reversal  is  occurring,  then  it  is  possible  that  the 
density  being  used  in  the  momentum  equation  will  be  varying  greatly  between  iterations.  This 
can  result  in  an  oscillation.  By  allowing  the  density  term  to  vary  from  the  upwind  donor  density 
to  a  fixed  donor  density,  large  changes  in  the  new  velocity  solution  are  reduced.  In  this  manner 
the  negative  impact  of  oscillating  flows  on  convergence  can  be  minimized.  This  donor  approach 
is  only  used  for  the  momentum  equation;  the  mass  and  energy  equations  always  use  the  upstream 
density.  Thus,  while  this  approach  may  result  in  errors  in  the  velocity  prediction,  it  will  not 
impact  the  overall  mass  and  energy  conservation. 

Along  with  the  conservation  equations,  an  equation  of  state  is  also  required.  This  is  given  by 


PV.  = 

X  1  T  1 


R  Ej 

Cp^(Ti) 


(6) 


In  the  above  equation,  Cp  is  the  specific  heat  of  the  compartment  at  its  current  tempera.ture, 
mwair  is  the  molecular  weight  of  air,  and  R  is  the  real  gas  constant. 


2.3  Conservation  Equations  in  Difference  Form 

Equations  1,3,  and  4  are  expressed  in  finite  difference  form  for  solution.  The  mass  and 
energy  equations  are  expressed  using  an  end-of-time  step  explicit  Euler  difference  scheme.  The 
momentum  equation  is  expressed  in  a  semi-implicit  form.  The  result  is  that  mass  and  energy 
will  be  conserved  absolutely;  however,  momentum  will  not  be  absolutely  conserved. 


Mass  Equation: 
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(7) 
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In  the  above  equation  n  indicates  the  current  time  step  number  and  n-1  the  prior  time  step 
withAt"  =  t"-t"'*. 

Energy  Equation: 


E"  = 


A 


2;a,v;F;AiP;h‘+E: 

V  j  J 


At"  +E"-' 


(8) 


Momentum  Equation: 


v”  =  v”  ^  + 

J  j  ^ 


At" 

PjLj 


(Pi"-p;+APj  + 


(pgAz)J-')- 


KjAt" 

~2L~ 


-  V 


(9) 


In  the  above  equation  n-  indicates  the  previous  iteration,  n+  indicates  a  next  time  step  gases, 
and  ~  indicates  a  linearized  prediction.  In  the  friction  term,  the  velocity  at  the  next  time  step  is 
approximated  with  a  tangent/secant  approach.  The  n+  value  is  the  velocity  from  the  prior 
iteration  if  the  velocity  direction  did  not  change  during  iterations,  and  it  is  zero  if  the  velocity 
changed  during  iterations. 


The  ~  values  for  the  pressure  indicate  an  approximation  of  the  next  time  step  pressure  based 
on  the  next  time  step  velocities.  This  approximation  accormts  for  the  change  in  pressure 
resulting  from  energy  inputs  into  the  control  volume  from  either  sources  and  sinks  such  as  fires 
and  wall  heat  transfer,  or  from  inputs  and  outputs  due  to  mass  flow  in  flow  connections. 


pn  _  pn-1 


dP"^  dE" 

+  ^- - !-At" 

dE"  dt 


(10) 
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(11) 

(12) 
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.  pn-1  Pf 
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1  n+ 

Si 


dE" 
^At" 


niw,,ViCp"^  dt 


(13) 


where  I  is  any  control  volume. 

Ideally,  we  want  to  use  the  new  specific  heat  in  determining  the  pressures  change.  However, 
that  involves  knowing  all  the  mass  flows  in  and  out  of  the  volume.  To  avoid  making  the  solution 
matrix  overly  large,  the  new  compartment  specific  heat  can  be  estimated  based  on  the  prior  time 
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step  velocities.  As  the  velocity  solution  converges,  the  specific  heat  solution  will  converge  as 
well.  Since  the  specific  heat  of  air  is  weakly  varying  with  temperature,  as  long  as  the  velocity 
and  resulting  temperature  change  in  a  compartment  are  limited  by  the  time  step  size,  the  lag  in 
specific  heat  convergence  will  have  a  minimal  impact  on  the  solution.  Since  the  linearized 
pressure  is  being  computed  based  on  an  energy  derivative,  in  order  to  conserve  energy  the  old 
pressure  must  be  multiplied  by  the  ratio  of  old  to  predicted  specific  heats  in  Equation  13.  Once 
the  extrapolated  pressures  are  computed,  the  momentum  equation  becomes: 
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In  the  case  of  a  duct  node  the  critical  equation  is  the  conservation  of  mass  since  mass  is  not 
stored  in  ducts.  Therefore,  for  a  duct  node.  Equation  7  is  solved  if  it  is  an  internal  node  that 
connects  multiple  ducts.  For  an  external  duct  node,  one  that  connects  a  duct  to  a  compartment, 
the  node  pressure  is  the  same  as  the  compartment  pressure;  therefore,  external  node  pressure  is 

not  a  solution  variable.  Thus  when  Equation  9  is  used  for  a  duct,  the  P”  ’s  will  be  replaced  by  a 
P"  for  each  internal  node. 


2.4  Hydraulic  Submodels 

2.4.1  Energy  Somce  Term 

The  E"  term  in  Equation  8  is  composed  of  a  number  of  energy  sources  and  sinks.  This  term 
is  given  by: 

E;  =  qf  (l  -  Xf  )-  qf.pyro  +  mf.pyrohf.pyro  “  qi,supp  +  mi,supphi.supp  +  qi,gas  "  S  (15) 

w 

where  q^  (l  -  Xf )  is  the  convective  heat  release,  qf  pyro  is  the  heat  of  pyrolysis  of  the 
fuel,  mj.py,.„h,.py^g  is  the  enthalpy  of  the  pyrolyzed  fuel,  qj  ^^p  is  the  energy  removal  due  to 
suppression  (e.g.  water  evaporation),  riijj^pphj  ^^p  is  the  enthalpy  of  the  suppression  agent  as  it  is 
added  to  the  control  volume,  is  the  radiant  energy  absorbed  or  loss  by  the  gas,  andq,-  ^  ^  is 
the  convective  heat  loss  to  a  surface. 
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2.4.2  Buoyancy  Head 

The  (pgAZ)  represents  the  buoyancy  head  through  a  junction  or  in  a  duct.  It  represents  the 
head  or  force  that  results  from  any  density  gradient  between  the  two  endpoints  of  the  junction  or 
duct.  For  a  junction  the  buoyancy  head  is  given  by: 


(pgAZ)j  =  j^(p, +^AZjj  +  Pi[^Zi  +^Az^ 


f 


■Pk 


1  . 

Zk+-Az, 


(16) 


For  a  duct  the  buoyancy  head  is  given  by: 
(pgAZ),  =p,(z^-z„,)g 


Pi 


If  either  node  nl  or  n2  is  a  terminal,  then  the  compartment  density  gradient  term, 

^  1  "l 

Z;  + — AZj  ,  must  be  included  in  the  manner  of  Equation  1 6. 

2  j 


(17) 


2.4.3  Loss  Terms 


Pressure  losses  occur  injunctions  and  in  ducts  and  duct  nodes.  Junction  pressure  losses  are 
estimated  as  a  discharge  coefficient  whose  value  is  a  function  of  the  type  of  opening.  These 
losses  are  represented  by  the  K  in  the  momentum  equation.  It  is  assumed  that  direction  through  a 
jimction  does  not  impact  the  flow  loss.  Thus,  only  one  K  is  needed.  Ducts  are  more  complex  as 
they  tend  to  contain  fittings  whose  loss  is  direction  dependent  such  as  an  expansion  joint  or  a  tee. 
Ducts  also  have  a  wetted  length  that  tends  to  be  much  larger  than  the  depth  of  a  typical  door  like 
opening.  Thus,  there  is  the  need  to  compute  a  wall  fnction  loss.  The  basic  equation  for 
determining  the  net  pressure  loss  in  a  duct  is  given  by  (Crane,  1988): 


AP. 


loss,d 


=  T? 


4fL„ 


(18) 


where  f  is  a  wall  fnction  factor,  s  is  a  duct  segment  of  constant  diameter,  and  i  is  a  fitting.  The 
fnction  factor  is  determined  by  the  duct  segment’s  Re3molds  number  and  wall  roughness,  s,  as 
determined  by  the  Colebrook  Equation.  As  this  equation  requires  an  iterative  solution,  a  direct 
solution  approximation  is  used  (George,  et  al.,  2000a): 
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(19) 


In  FSSIM,  ducts  have  unidirectional  flow.  Thus,  as  long  as  a  length  of  duct  is  not 
interrupted  by  either  terminating  in  a  compartment  or  in  a  fitting  with  multiple  other  ducts  (e.g. 
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where  more  than  one  exit  flow  path  exists),  this  can  be  considered  one  duct  for  the  purpose  of  the 
hydraulic  calculation.  However,  this  means  that  one  duct  could  have  multiple  cross  sectional 
areas.  To  avoid  having  to  store  the  duct  segment  information,  an  effective  loss  can  be 
determined.  The  loss  is  a  function  of  velocity,  and  in  the  case  of  the  wall  friction,  the  diameter 
as  well.  Since  a  duct  has  unidirectional  flow  and  no  mass  storage,  the  volumetric  flow  rate 
through  any  cross  section  of  a  duct  is  constant.  Thus,  all  losses  can  be  expressed  in  terms  of  a 
single  diameter.  The  following  derivation  is  shown  to  demonstrate  how  this  is  accomplished. 
FSSIM  does  not  perform  this  calculation  internally,  but  rather  the  user  would  need  to  do  this 
when  defining  the  duct  for  input. 


=  K 


=  K. 


eq 


eq 


feq^s  2 

- 2 - 

d 

eq 


eq 


fs^eq  Vs  ^  fs^eq 
ds  d,  A3 


(20) 

(21) 


The  friction  factors  would  still  require  computation  at  each  time  step,  as  they  are  dependent 
on  a  complex  relationship  with  the  Reynolds  number.  However,  if  fully  turbulent  flow  is 
assumed  in  the  duct,  then  f  depends  only  on  the  wall  roughness  and  duct  diameter.  With  this 
assumption: 


^^loss,d  ~ 


_Pd 


K  V 


4fL. 


eq 


(22) 


In  addition  to  losses  along  a  contiguous  length  of  duct,  there  are  losses  associated  with  the 
joining  of  multiple  duct  lengths  such  as  a  tee  or  manifold.  If  there  are  n  ducts  connected  to  node, 
then  there  will  be  n(n-l)  possible  combinations  of  inlet  and  outlet  flow.  For  each  incoming  flow 
direction  there  is  a  loss  associated  with  that  incoming  flow  exiting  the  other  flow  directions 
given  n(n-l)  form  losses.  Computing  the  added  loss  in  a  ductnode  is  done  before  each  iteration 
based  on  the  prior  iteration’s  velocity  and  density  solution.  First,  each  inlet  duct  is  assigned  a 
weighting  factor  based  on  its  incoming  flow  fraction.  Each  outlet  duct  is  then  assigned  a  loss 
based  on  the  weighted  sums  of  the  incoming  flow  loss  times  the  appropriate  form  loss  factor. 

The  calculated  equivalent  loss,  Kd,  is  used  for  Kj  in  Equation  9. 


2.4.4  Bidirectional  Flow 


Two  submodels  for  bidirectional  flow  are  used  in  FSSIM.  The  first  model  is  for  horizontal 
flow  through  vertical  openings  (e.g.  doors).  The  second  model  is  for  vertical  flow  through 
horizontal  openings  (e.g.  hatches).  The  horizontal  flow  model  uses  the  neutral  plane  concept  and 
the  assumption  of  linear  pressure  gradients  within  compartments  (Karlsson  and  Quintiere,  2000) 
to  determine  at  the  beginning  of  each  time  step  how  a  junction  is  to  be  partitioned.  The  vertical 
flow  model  is  based  on  the  work  performed  by  Cooper  (1989). 
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2. 4. 4. 1  Bidirectional  Horizontal  Flow  Submodel 

Horizontal  flow  in  an  opening  between  compartments  is  driven  by  differences  in  relative 
pressure  over  the  height  of  the  opening.  If  each  compartment’s  pressure  is  known  at  its  floor  and 
each  compartment  is  at  a  rmiform  density,  then  the  pressure  as  a  function  of  height  for  each 
compartment  can  be  given  as: 

Pi(z)  =  Pi(zo)-Pig(z-Zo)  (23) 

The  differential  pressure  in  the  opening  is  given  as: 

APj.(z)  =  Pi  (z-e )-  Pk  (zk,o )+  Pkg(z  -  z,^o )-  Pig(z  -  z.  0 )  (24) 


If  the  differential  pressure  is  positive,  then  flow  will  be  from  compartment  i  to  k  at  that 
elevation,  and  the  flow  will  be  in  the  reverse  direction  if  the  differential  pressme  is  negative.  If 
the  sign  of  the  differential  pressure  changes  over  the  height  of  an  opening,  then  there  will  be 
bidirectional  flow.  The  location  of  the  neutral  plane,  the  location  where  the  pressure  differential 
is  zero,  is  given  by: 

„  _  Pk(zk,o)-Pi(zi,o)+g(PkZk.o  -PiZi.o) 

Fra 

Each  horizontal  flow  junction  is  represented  by  the  model  as  two  parallel  junctions.  By 
applying  Equation  25  at  the  beginning  of  each  time  step  the  flow  acceleration  through  the 
junction  can  be  partitioned  for  bidirectional  flow.  To  preserve  the  prior  time  step  mass  flow 
through  the  junction  the  velocities  at  the  beginning  of  a  times  tep  need  to  be  adjusted  by  the 
computed  change  in  flow  area.  To  avoid  numerical  stability  problems,  the  junction  area  change 
is  relaxed  and  neither  parallel  junction  is  allowed  to  be  less  than  10  %  of  the  flow  area. 

2. 4. 4. 2  Bidirectional  Vertical  Flow  Submodel 

A  strict  application  of  Equation  4  would  result  in  all  horizontal  openings  having 
unidirectional  flow.  Since  each  compartment  is  represented  by  one  set  of  quantities  and  the 
opening  has  no  elevation  change,  as  it  is  horizontal,  the  pressure  change  across  the  opening  is 
uniform.  In  reality,  however,  imstable  flow  conditions  may  exist  at  these  openings  if  cold  dense 
air  is  above  warmer  less  dense  air.  If  denser  gas  is  located  in  the  upper  compartment,  depending 
upon  the  vent  conditions  there  may  be  flow  in  both  directions.  Bidirectional  flow  will  occur 
through  a  horizontal  opening  unless  the  pressure  difference  exceeds  a  flooding  criterion.  This 
criteria,  where  Cs  is  a  shape  factor  equal  to  0.754  for  a  circular  opening,  is  (Cooper,  1989): 
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(26) 
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Cooper’s  correlations  for  vertical  flow  involve  a  three-step  process.  First,  the  flowrate  is 
determined  as  normal  (e.g.  assuming  unidirectional  flow).  Second,  the  flooding  limit  is 
determined.  Third,  if  the  flooding  limit  has  not  been  exceeded  an  exchange  flowrate  is 
calculated  along  with  a  new  flowrate  in  the  opposite  direction  such  that  the  sum  of  the  flowrates 
is  equal  to  the  first  step.  The  exchange  flowrate  is  given  as: 
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V-  =  V- 

J»ex  j, ex  max 


1- 


Iap, 
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AP. 


r  V- 

ex  max 


(27) 


The  Vj,ex  in  the  above  equation  represents  a  steady-state  velocity.  Equation  4  states  that 
steady  state  will  be  reached  when  the  pressure  across  an  opening  equals  the  loss  through  the 
opening.  Thus,  Equations  4  and  27  can  be  combined  to  yield  a  finite  differenced,  time  derivative 
form  of  the  exchange  flow: 


V- 

2L, 


v” 

J,ex 


v"'" 

J.ex 


v"  =  v"  '  + 

j,ex  j,ex 


At”  1^ 

- Kj 

Lj  2  ' 


0.055 


4F;Aj  j8gAgjF;Aj 


71  y7i(pi+pj 


1-1 


AR 


AP, 


j,  flood 


+  • 


K,At" 


J) 


2L, 


J.ex 


v” 

J.ex 


(28) 


The  aj,ex  is  a  direction  function  that  will  allow  for  flow  acceleration  along  the  proper 
direction  for  the  exchange  flow.  For  flow  in  the  primary  direction  the  vj  terms  in  Equation  14  are 
replaced  by  vj  +  Vj,ex- 

2.4.5  Specific  Heat 

The  equation  of  state  used  is  the  ideal  gas  law  expressed  in  energy.  Equation  6.  This 
equation  relates  the  enthalpy  of  the  gases  in  a  compartment  to  the  compartment  temperature  by 
means  of  the  specific  heat: 

E."=M”T”Cp(Ti")  (29) 

The  specific  heat  is  a  non-hnear  function  of  temperature.  Thus,  determining  the  value  of  the 
specific  heat  required  for  a  particular  compartment  enthalpy,  which  is  a  conserved  quantity,  and 
pressure  state  requires  simultaneously  determining  the  compartment  temperature.  This  can  be 
done  in  an  iterative  fashion.  A  guess  at  the  compartment  temperature  is  made  using  the  prior 
value  for  Cp  and  the  current  enthalpy.  This  temperature  is  used  to  re-evaluate  Cp  using  the 
Newton-Raphson  method  (Hildebrand,  1976),  and  the  process  is  repeated  until  Cp  and 
temperature  converge.  This  same  process  can  be  used  for  determining  the  effective  temperature 
of  a  duct  node  given  its  pressure  and  incoming  flows. 
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2.5  Combustion 


2.5.1  Ignition 

As  part  of  the  FSSIM  input,  users  can  define  one  or  more  fires  to  start  at  explicit  times. 
User’s  may  also  choose  to  allow  fire  to  spread.  Ignition  of  additional  fires  is  determined  at  the 
beginning  of  each  time  step.  Each  compartment  can  have  a  usetype  designated  for  it,  which 
denotes  a  fuel  loading  and  fuel  classification.  Separate  temperature  ignition  criteria  can  be 
defined  for  surfaces,  temperature  of  incoming  vent  flows,  and  compartment  temperature  (Back, 
et  al.,  2003).  Overhead  surfaces  can  be  given  a  different  ignition  temperature  from  non-overhead 
surfaces. 

2.5.2  Pyrolysis 

Pyrolysis  is  determined  by  one  of  three  methods:  the  fire  has  a  constant  pyrolysis  rate,  a  t^ 
pyrolysis  rate,  or  a  user-defined  pyrolysis  rate.  Growth  in  pyrolysis  can  be  limited  by  specifying 
a  maximum  pyrolysis  rate  in  kg/m^-s.  All  fires  can  be  given  an  end  time  in  either  absolute  time 
or  fuel  loading.  The  calculated  pyrolysis  rate  can  be  reduced  be  various  mechanisms.  If  the  fire 
is  being  suppressed  by  an  agent,  the  pyrolysis  will  be  reduced  by  a  suppression  factor.  If  the  fire 
has  become  oxygen  limited,  then  the  pyrolysis  rate  is  determined  by  a  linear  function  of 
temperature.  The  maximum  pyrolysis  rate  after  oxygen-limiting  conditions  are  reached  is  set  at 
the  point  where  the  fire  became  oxygen  limited;  the  compartment  temperature  for  that  point  is 
also  stored.  The  oxygen-limited  pyrolysis  rate  is  then  calculated  by  the  ratio  of  the  current 
temperature  to  the  stored  temperature  not  to  exceed  the  maximum  rate. 

2.5.3  Combustion  Heat  Release  +  Species  Production 

Combustion  of  pyrolyzed  fuel  is  calculated  on  the  basis  of  the  available  oxygen  in  the 
compartment  where  the  pyrolysis  is  occurring.  Currently,  there  is  no  burning  of  xmbumt  fuel 
exiting  one  oxygen-depleted  compartment  into  a  second  oxygen-rich  compartment.  If  there  is 
sufficient  oxygen  in  the  compartment  above  a  user  specified  lower  limit,  all  the  pyrolyzed  fuel 
will  combust.  If  there  is  insufficient  oxygen  in  the  compartment,  a  more  detailed  estimate  of  the 
available  oxygen  is  made  which  includes  a  prediction  of  the  net  inflow  of  oxygen  based  on  prior 
time  step  velocities.  The  actual  heat  release  rate  is  adjusted  to  use  the  calculated  amount  of 
available  oxygen.  Note  that  a  fire  in  an  oxygen  depleted  compartment  can  bum  at  a  rate  equal  to 
the  pyrolysis  if  there  is  a  sufficient  flowrate  of  oxygen  into  the  compartment. 

Species  are  generated  based  on  user-provided  yields  for  each  fuel  being  burned.  These 
)delds  represent  the  mass  of  combustion  products  formed  for  a  imit  mass  of  fuel  burned.  Note 
that  the  consumption  of  oxygen  can  be  expressed  as  a  negative  yield. 

2.6  Overall  Solution  Scheme 

An  iterative  solution  scheme  is  used  to  solve  the  hydraulic  equations.  This  scheme  consists 
of  time  step  initialization,  an  outer  loop,  an  inner  loop,  and  a  time  step  wrap-up.  The  outer  loop 
goal  is  convergence  of  the  compartment  pressure  solution  and  prevention  of  large  changes  in 
thermophysical  conditions  over  a  time  step.  The  inner  loop  goal  is  convergence  of  the  junction 
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and  duct  velocity  solution  along  with  conservation  of  mass  at  duct  nodes.  The  flow  chart  in 
Figure  2  shows  the  iteration  logic. 


Figure  2:  Iteration  Flowchart 


2.6.1  Time  step  Initialization 

Before  entering  the  outer  loop,  the  time  step  initial  conditions  are  established.  Output  is 
dumped  if  a  user-defined  output  interval  has  elapsed.  All  control  functions  are  evaluated  using  a 
recursive  routine;  the  recursion  allows  for  control  functions  that  are  functions  of  other  control 
functions.  The  available  flow  area  of  jimctions  and  ducts  are  determined.  Finally,  compartments 
are  checked  for  ignition  of  a  new  fire,  and  any  existing  fire’s  pyrolysis  rate  is  determined. 

2.6.2  Outer  Loop 

A  time  step  starts  outside  the  outer  loop.  End  of  time  step  Cp  values  for  each  compartment 
are  predicted  from  the  prior  time  step  solution.  The  inner  loop  is  then  executed.  Upon  exiting 

the  inner  loop,  the  P;  values  used  in  Equation  14  are  compared  to  solution  values,  P” .  If  the 
difference  in  value  is  below  a  convergence  criteria,  the  time  step  is  advanced.  This  step  ensures 
that  linearized  pressure  relations  for  mass  flows  used  in  the  momentum  equation  were  reasonable 

at  the  end  of  the  time  step.  If  not,  the  outer  loop  is  cycled  and  the  Pj  values  will  be 
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approximated  using  the  prior  iteration’s  solution.  If  the  P;  values  are  converged,  then  the 
velocities,  temperatures,  and  pressures  are  checked  for  their  relative  change  and  a  Courant 
criterion  is  checked.  If  any  of  the  associated  limits,  to  be  discussed  later,  are  exceeded  the  time 
step  is  reduced  and  the  outer  loop  is  cycled.  Lastly,  the  new  time  step  size  is  determined  by 
determining  the  variable  that  was  closest  to  exceeding  its  convergence  criteria.  The  new  time 
step  is  given  a  multipher  and  the  current  time  step.  The  multiplier  is  the  ratio  of  the  closest 
convergence  to  its  associated  criteria  not  to  exceed  a  user-defined  maximum  time  step.  If  either 
too  many  outer  loops  iterations  or  too  small  of  a  time  step  occurs,  then  the  time  step  will  fail  and 
the  calculation  will  abort. 

2.6.3  Inner  Loop 

The  inner  loop  constructs  and  solves  the  matrix  for  the  duct  and  junction  momentum 
equations  and  mass  conservation  at  internal  duct  nodes.  However,  since  the  new  velocities  are  a 
fimction  of  the  energy  source  terms  in  each  compartment,  the  inner  loop  also  updates  the  somce 
terms  for  each  iteration.  The  first  step  in  the  inner  loop  is  setting  the  v"^  values  by  examining 
ducts  and  junctions  for  flow  reversal.  Second,  donor  values  are  determined  using  Equation  5. 

The  prior  iteration  velocities  and  source  terms  are  used  to  generate  the  Cp""^  values  in  Equation 

13.  The  compartment  heat  release  rate  is  then  determined  along  with  the  effects  of  suppression 
systems.  Convective  and  radiative  heat  transfer  are  then  determined  for  each  compartment.  The 
solution  matrix  is  constructed,  zero-value  columns  and  rows  are  eliminated  (resulting  firom  ducts 
and  jimctions  that  have  dampers  or  doors  closed),  and  the  matrix  solved.  Equations  7  and  8  are 
used  to  update  compartments  and  duct  nodes.  The  velocities  and  duct  node  mass  conservation 
are  checked  for  convergence,  criteria  to  be  discussed  later.  If  the  criteria  are  not  met,  the  inner 
loop  is  cycled.  If  the  number  of  inner  iterations  exceeds  a  user  specified  criteria,  the  inner  loop 
is  exited,  the  time  step  reduced,  and  the  outer  loop  is  cycled. 

If  the  inner  loop  is  cycled  a  partial  update  is  performed  of  the  energy  source  terms.  Heat 
release  rates  are  updated  in  ventilation-controlled  compartments  using  the  new  velocity  solution; 
this  is  not  needed  in  non-ventilation  controlled  compartments  as  their  heat  release  rate  is  not  a 
function  of  incoming  mass  flows.  The  radiation  solution  is  updated  in  compartments  containing 
a  fire  using  the  newly  predicted  compartment  and  surface  temperatures  resulting  from  the  current 
velocity  solution.  Non-burning  compartments  are  not  updated  in  order  to  reduce  computational 
time;  their  radiation  source  term  is  not  likely  to  vary  significantly  from  the  first  iteration. 

Lastly,  surface  heat  fluxes  are  recomputed. 

2.6.4  Time  step  Wrap-up 

After  successful  completion  of  the  outer  loop,  the  time  step  is  advanced.  A  final  update  is 
performed  on  the  surface  temperatures.  All  other  variables  are  also  updated  by  setting  all  prior 
time  step  values  to  the  newly  calculated  values.  Lastly,  the  simulation  clock  is  advanced. 
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2.6.5  Convergence  Criteria 

2. 6. 5. 1  Inner  Convergence  Criteria 

There  are  three  sets  of  convergence  criteria  for  the  inner  iterations.  The  first  criterion  is  for 
the  mass  conservation  for  a  duct  node.  The  second  set  is  for  the  velocity  in  ducts.  The  third  set 
of  criteria  is  for  velocity  injunctions.  Mass  conservation  is  determined  by  calculating  the  net 
mass  flow  at  a  duct  node  and  the  sum  of  the  absolute  value  of  the  mass  flows  at  a  duct  node.  If 
conservation  is  met,  the  net  mass  flow  should  be  zero;  however,  round  off  error  in  the  floating 
point  math  may  prevent  this.  Instead,  the  net  mass  flow  is  compared  to  the  total  mass  flow,  and 
if  the  error  is  greater  than  0.001  %,  the  inner  loop  is  cycled.  The  duct  and  junction  velocity 
criteria  have  similar  logic  but  can  be  set  with  different  values  for  convergence.  For  each  duct 
and  junction,  if  the  flow  area  was  nonzero,  and  the  velocity  magnitude  of  the  prior  time  step  and 
current  solution  are  both  above  a  cutoff,  currently  1  pm/s,  convergence  is  checked.  Without  the 
cutoff  junctions  and  ducts  with  essentially  zero  flow  that  have  no  impact  on  the  overall  mass 
exchange  could  prevent  convergence.  Two  criteria  are  then  checked.  The  first  is  flow  reversal; 
if  the  flow  in  a  duct  or  junction  changes  direction,  the  iimer  iteration  is  repeated.  The  second  is 
the  relative  change  in  velocity.  If  the  relative  change  in  velocity  between  iterations  is  greater 
than  5%,  the  inner  iteration  is  cycled. 

2. 6. 5. 2  Outer  Convergence  Criteria 

There  are  two  sets  of  convergence  criteria  for  the  outer  loop.  The  first  compares  the 
Pj  values  used  in  Equation  14  to  the  solution  values,  P” .  If  they  differ  by  more  than  1  %,  the 
outer  iteration  is  repeated  without  a  change  in  time  step.  This  criterion  ensures  that  the 
thermophysical  properties  assumed  during  the  solution  reflect  the  end  of  time  step  values  for 
those  properties.  The  second  set  of  criteria  examines  the  relative  change  in  temperature,  mass 
flows,  and  pressures.  Specifically,  if  temperature  in  any  compartment  changes  by  more  than 
10  %,  pressure  in  any  compartment  changes  by  more  than  1  %,  or  the  net  mass  flow  for  a 
compartment  exceeds  50  %  of  the  end  of  beginning  of  time  step  mass,  then  the  outer  loop  is 
cycled  with  a  reduced  time  step.  The  time  step  reduction  is  the  larger  of  50  %  or  the  largest 
percentage  by  which  any  convergence  criterion  was  exceeded. 
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3.0  HEAT  TRANSFER  SOLVERS 
3.1  Convection 


Convection  heat  transfer  is  calculated  for  each  surface  for  each  time  step  using  the  beginning 
of  time  step  surface  and  compartment  temperatures.  Depending  on  the  surface  orientation, 
horizontal  or  vertical,  one  of  two  natural  convection  heat  transfer  correlations  are  used  to 
determine  the  convective  heat  transfer.  The  correlations  used  are  (Holman,  1990): 


Horizontal:  heated  facing  downward  or  cooled  facing  upward 


C  =  0.59 


1/4 


(30) 


Horizontal:  heated  facing  upward  or  cooled  facing  downward 

q";"  =  1.52  iTi"  -  T:  r  (Ti"  -  )  (31) 


Vertical 
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where  i  indicates  a  compartment,  L  is  a  characteristic  length  of  the  surface  (taken  as  the  square 
root  of  the  surface  area),  and  w  indicates  a  surface. 

3.2  Radiation 

Radiation  heat  transfer  is  computed  on  a  compartment-by-compartment  basis  using  the 
beginning  of  time  step  surface  temperatures,  compartment  temperature,  compartment  gas 
composition,  and  compartment  heat  release  rate.  The  heat  transfer  is  calculated  using  a 
modified  gray-gas,  n-surface  net  radiation  method  (Forney,  1991). 

3.2.1  Net  Radiation  Equation  Derivation 


The  radiation  heat  transfer  equation  for  a  smface  is: 
A^oEwT^J +(l-8„)q"  =  q”  +  A„Aq'; 


(33) 


with 

A  Aq"  =q°“‘-q” 

^  *-w  Mw  'tw 


(34) 
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C  =A^as„T^+(l-8„)q" 


*iw  ~  ^4w2^w2-w^w2-w  ■*■  ’ 


w2 


(35) 

(36) 


where  Fw2-w  is  the  view  factor  from  surface  w2  to  w.  In  the  absence  of  attenuation  or  scatter  this 
is  the  fraction  of  energy  leaving  w2  that  arrives  at  w.  t  is  the  transmission  factor  from  surface 
w2  to  w,  the  fraction  of  energy  transmitted  over  the  distance  between  surface  w2  and  w.  Cw  is  a 
source  term  representing  such  as  radiant  emissions  from  a  fire  or  hot  gasses  in  the  compartment. 
Aq  represents  the  net  radiant  heat  transfer,  and  8  is  the  surface  emissivity. 

If  Equations  33  through  36  are  combined  in  conjunction  with  the  view  factor  symmetry 
relationship,  AwFw-w2  =  Aw2Fw2-w,  an  equation  for  the  net  heat  transfer  can  be  derived. 


w2“W 


w2 


^w2 


'J^w2-w  =  -2]^T^Fw-w2'^w-w2 


w2 


(37) 


Since  it  is  numerically  possible  for  the  emissivity  of  a  surface  to  be  zero,  the  above  equation 
is  modified  by  the  substitution 


Aq;  =  8„Aq; 


(38) 


Aq;  -li(l“^w2)Aq;2Fw2-w'Cw2-w 

w2 


w2 


(39) 


This  results  in  a  dense  matrix,  the  only  zero  terms  are  for  surfaces  along  the  same  wall, 
which  can  be  solved  by  using  a  number  of  numerical  techniques.  This  matrix  will,  however, 
always  be  diagonally  dominant  since  0  <  8  <  1.  The  transmission  factor  is  calculated  by 
assuming  a  gray  gas  and  it  is  given  by 


<7-  ^  £i““Lw2-w 

‘'W2-W  ”  ^ 


(40) 


where  a  is  an  absorption  coefficient  which  is  calculated  based  on  the  amounts  of  CO2,  H2O,  and 
soot  in  the  compartment  along  with  an  average  pathlength.  L  represents  an  average  path  length 
between  two  radiatively  communicating  surfaces.  Note  that  this  is  not  the  same  pathlength  used 
for  calculating  a.  a  is  calculated  using  the  ABSORB  subroutine  from  CFAST  v3.1.7  (Jones,  et 
al.,  2000).  This  routine  calculates  a  by  doing  a  table  lookup  of  absorptivities  of  the  CO2  and 
H2O  mass  fractions  and  accounting  for  band  overlap,  which  can  be  approximated  as  one  half  of 
the  CO2  absorptivity  (Tien,  Lee,  and  Stretton,  2002).  The  net  a  is  calculated  by  summing  the 
table  values  with  a  contribution  due  to  soot: 


•^gas  ~  ^^HjO  ^C02  2 


(41) 
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(42) 


a  =  1  4k  V  T  - 


where  ks  is  an  empirical  constant  (ks==l  195),  Vi,s  is  the  soot  volume  fraction  in  compartment  i,  T, 
is  the  gas  temperature,  and  U  is  the  beam  length  calculated  with  Equation  44. 

The  source  term,  Cw,  is  determined  by  assuming  a  gray  gas  and  that  the  fire  radiates  a  fixed 
fraction  of  its  heat  release  as  radiant  energy. 

qf  +  E "  '^w2-w  ,  (43) 

w2 


where  %  is  the  radiant  fraction  of  the  fire. 


3.2.2  Solution  of  the  Net  Radiation  Equation 

To  simplify  the  solution  of  equation  39  three  assumptions  are  made.  First,  each  compartment 
is  assumed  to  be  a  rectangular  parallelepiped.  Every  surface  in  each  compartment  is  then 
associated  with  one  of  the  six  surfaces  of  the  parallelepiped.  Second,  the  view  factor  of  any 
subsurface  for  a  wall  of  the  parallelepiped  is  determined  from  the  area  fraction  of  that  surface 
with  respect  to  the  parallelepiped  wall  containing  it.  Third,  the  fire  is  located  at  the  center  of  the 
parallelepiped  and  “paints”  the  walls  based  on  their  area  fraction  with  respect  to  the  total 
compartment  surface  area.  The  calculation  of  pathlengths  and  view  factors  is  greatly  simplified 
by  these  assumptions. 


Pathlengths  for  fire  radiation  to  the  surfaces  is  calculated  as  one  half  of  the  compartment 
dimensions.  Thus,  the  pathlength  from  the  fire  to  the  floor  or  ceiling  is  given  as  one  half  the 
height  of  the  compartment.  Pathlengths  for  surface-to-surface  communication  is  calculated  as 
the  parallelepiped  wall  center-to-wall  center  distance.  Thus,  the  floor-to-ceiling  pathlength  is  the 
height  of  the  compartment  and  the  ceiling-to-sidewall  pathlength  would  be  the  diagonal  line 
connecting  the  ceiling  center  to  the  side  center.  For  the  purpose  of  calculating  absorption 
coefficients  the  pathlength  is  calculated  as  shown  below  (Holman,  1990): 


Li  = 


3.6Vi 

Za. 


w 


(44) 


View  factors  from  the  fire  to  the  surfaces  are  calculated  solely  on  the  basis  of  the  area 
fraction  of  the  wall  with  respect  to  the  entire  compartment. 


Ff- 


(45) 


Surface  view  factors  are  calculated  as  either  parallel-plate  view  factors  or  right-angle-comer 
view  factors  based  on  the  parallelepiped  walls  to  which  the  surfaces  belong  to.  The  view  factor 
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of  each  surface  is  given  as  the  appropriate  wall  view  factor  times  the  area  fraction  of  the  surface 
with  respect  to  that  wall.  Parallel-plate  view  factors  are  given  as: 


F  = 


Tixy 
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where  x  and  y  are  the  ratios  of  each  side  to  the  plate  separation.  Comer  view  factors  are  given 
as: 
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where  z  is  the  length  of  the  common  side,  x  and  y  are  the  lengths  of  the  other  two  sides  of 
rectangles  1  and  2,  respectively,  X=  x/z,  and  Y=y/z  (Siegel,  et  al.,  2001) 

The  first  part  of  the  solution  is  to  consfruct  the  right  hand  side  vector  for  the  net  radiation 
matrix.  This  vector  consists  of  the  source  term  representing  radiative  inputs  from  the  fire  and  the 
gas  as  well  as  the  direct  emissions  from  the  surfaces;  reflected  emissions  are  part  of  the  net 
radiation  terms  on  the  left-hand  side. 
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where  side(w)  indicates  the  compartment  face  that  contains  surface  w. 
The  left-hand  side  is  consfructed  as  follows: 
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After  solving  the  matrix  the  net  radiation  is  extracted  from  the  modified  net  radiation.  In 
both  Equations  48  and  49  the  treatment  of  the  gas  and  surface  temperatures  depends  on  the 
current  inner  iteration.  For  the  first  iteration,  they  are  taken  from  the  previous  time  step’s 
solution.  For  subsequent  inner  iterations,  they  are  taken  as  the  average  of  the  previous  time  step 
and  the  previous  iteration.  At  high  temperatures,  small  changes  value  can  have  large  changes  in 
the  radiative  field.  By  taking  the  average  of  the  prior  time  step  and  the  prior  iteration,  this  effect 
can  be  minimized. 


Aq  =  Aq *  8  A 


(50) 


The  net  radiation  absorption  by  the  gas  is  determined  by  the  difference  in  the  total  net 
radiation  to  the  surfaces  and  the  radiative  output  of  the  fire.  Any  difference  between  these  two 
values  results  from  a  net  emission  or  a  net  absorption  by  the  gas  in  the  compartment. 

Aqi  =Xfqf+2]AqXA,  (51) 


3.3  Conduction 

3.3.1  ID  Conduction  Equation 

Conduction  heat  transfer  is  computed  for  each  surface  at  the  end  of  each  time  step  using  the 
beginning  of  time  step  computed  convective  and  radiative  heat  transfer  solutions  as  boundary 
conditions.  Heat  conduction  is  calculated  using  a  ID  solver  with  multiple  layers  of  materials  and 
temperature-dependent  specific  heats  and  conductivities.  The  general  equation  for  ID  heat 
transfer  in  Cartesian  coordinates  is  (Holman,  1990): 


— k— 
5x  dx 


+q 


(52) 


where  k  is  the  conductivity,  p  is  the  material  density,  q^is  a  volumetric  heat  generation  rate,  T  is 
the  material  temperature,  x  is  the  position  within  the  material,  and  Cp  is  the  material  specific  heat. 
In  general  k,  p,  and  c  are  all  functions  of  temperature.  However,  allowing  p  to  be  a  function  of 
temperature  means  the  solution  would  either  have  to  take  into  account  expansion  of  the  material 
which  adds  greatly  to  the  required  computation  resources  required,  or  conservation  of  mass 
would  not  be  maintained  which  is  also  undesirable.  At  the  boundaries  the  following  holds: 


(53) 


where  w  is  the  material  surface  and  q"  is  the  net  surface  heat  flux. 
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3.3.2  Solution  of  ID  Conduction  Equation 


To  solve  Equation  52  a  few  simplifications  are  made:  no  internal  heat  generation,  no  removal 
of  material  from  a  surface  (e.g.  as  would  represent  mass  loss  from  pyrolysis),  and  constant 
density  (though  there  is  temperature-dependent  conductivity  and  specific  heat).  The  equation  is 
discretized  with  central  differences  in  space  and  with  an  implicit  Crank-Nicholson  scheme  in 
time  (Strauss,  1992).  The  ID  surface  is  divided  into  N-1  cells  or  nodes  with  N  boundaries. 
Temperatures  are  stored  at  the  boundary  locations.  Thermophysical  properties  are  explicit  firom 
the  previous  time  step  solution.  The  general  form  of  the  spatial  discretization  is: 
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where  i  is  a  node  and  Ax  is  a  node  spacing.  Discretization  of  the  time  step  yields: 
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Boundary  conditions  are  implemented  by  assuming  one  of  two  boundary  conditions.  The 
first  time,  through  the  solver  the  boxmdary  condition  is  a  constant  radiative  heat  flux  and  a 
constant  convective  heat  transfer  coefficient  based  on  the  beginning  of  time  step  compartment 
temperature  and  radiative  conditions  from  Equations  30-32  and  50.  The  radiative  flux  is 
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corrected  for  the  new  surface  temperature  to  avoid  instabilities  at  high  surface  temperatures. 
This  correction  is  performed  by  linearizing  the  radiative  heat  transfer  change  as  follows: 
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If  the  time  step  decreases  during  an  iteration,  the  heat  transfer  solution  is  repeated  with  a 
constant  heat  flux  condition  equal  to  the  sum  of  the  constant  radiative  flux  and  the  time-average 
convective  flux  determined  from  the  first  pass  through  the  solver.  At  the  i=l  and  i=N  surfaces, 
the  general  discretized  equation  for  the  first  pass  through  the  solver  becomes: 
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Equations  56,  58,  and  59  result  in  a  tridiagonal  matrix.  The  matrix  is  solved  for  the  current 
time  step  using  a  tridiagonal  solver  and  the  relative  difference  in  the  new  solution  compared  to 
the  previous  solution  is  calculated.  The  node-spacing  values  are  selected  automatically.  The 
terminal  nodes  for  each  layer  of  a  surface  are  set  such  that  their  Biot  number  is  at  least  0.1  with 
the  heat  transfer  coefficient  assumed  to  be  100  W/m^-K.  Node  sizes  progressing  into  a  layer  are 
obtained  by  successive  doubling  of  the  outer  node  size  (for  example:  0.1  m,  0.2  m,  0.4  m,  0.2  m, 
0.1  m).  This  technique  minimizes  the  number  of  nodes  while  maintaining  a  reasonable  accuracy 
in  the  time-dependent  solution  (George,  et  al.,  2000b).  If  the  temperature  change  in  any  one  cell 
exceeds  a  specified  tolerance,  currently  0.1%,  the  solution  time  step  is  divided  by  an  integer 
factor  representing  the  factor  by  which  the  tolerance  was  exceeded.  The  matrix  is  now  solved 
within  a  loop,  with  the  thermophysical  parameters  calculated  at  the  start  of  each  sub-time  step. 

After  updating  a  surface,  a  correction  is  made  to  the  radiation  heat  transfer  to  the  gas.  The 
radiant  heat  transfer  boundary  conditions  in  Equations  58  and  59  account  for  the  change  in 
surface  temperature  over  a  time  step. 
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If  the  Biot  number  calculation  for  a  single-layer  surface  results  in  only  one  node  for  that 
surface,  the  above  method  is  not  used.  Instead,  a  lumped  parameter  analysis  is  performed  where: 

T"  =  T"-'  +  (q""'  +  q^"-' )  (60) 

Axpc 

4.0  MISCELLANEOUS 

4.1  Control  Functions 

FSSIM  has  a  limited  capability  to  model  control  systems.  Two  types  of  control  functions  are 
available.  The  first  type  is  a  binary  (on/off)  control  function  which  will  change  state  from  true  to 
false,  or  vice  versa,  when  a  user-specified  criterion,  such  as  time  or  a  compartment  temperature, 
is  reached.  This  function  can  be  defined  as  a  trip  function.  As  a  trip,  the  variable  will  only 
change  state  once,  so  if,  for  example,  a  temperature  criterion  is  reached,  the  function  will  not 
change  state  if  the  temperature  decreases  below  the  setpoint.  The  second  type  of  function  is  an 
AND  gate  which  will  change  the  state  of  a  binary  function  when  two  specified  conditions  are 
true.  This  function  type  can  be  used  to  specify  equipment  operation  based  on  some  physical 
criteria  along  with  an  equipment  state  variable  to  indicate  its  ability  to  function. 

4.2  Fire  Detection  Models 

A  number  of  detection  devices  are  modeled  by  FSSIM.  All  of  the  detection  devices,  with  the 
exception  of  UV  or  IR  detectors,  can  be  located  in  either  a  compartment  or  a  duct. 

4.2. 1  lonization/Photoelectric  Detector 

Ionization  and  photoelectric  detectors  will  alarm  based  on  a  user  provided  optical  depth;  the 
default  is  0.0581  m  '.  The  optical  depth  is  calculated  as  (Mulholland  and  Croarkin,  2000): 


OD  =  a3p,,  (61) 

where  ps  is  the  soot  density  and  as  is  the  mass  extinction  coefficient  of  soot,  8700  m  /kg. 

4.2.2  Gas  Detector 

A  gas  detector  will  alarm  when  the  gas  it  measures  reaches  a  user-defined  mass  fraction. 

4.2.3  Fixed  Temperature  Detector 

A  fixed  temperature  detector  will  alarm  when  the  temperature  it  measures  reaches  a 
user-defined  threshold. 
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4.2.4  RTI  Detector 


The  RTI  detector  type  simulates  the  response  of  a  thermal  activation  link  as  might  be  present 
on  a  sprinkler  head.  The  temperature  change  in  the  detector  is  given  by: 


where  Tj  or  j  is  the  temperature  of  the  compartment  or  duct  in  which  the  detector  is  located  and 
Vfor  j  is  either  a  ceiling  jet  velocity  based  on  the  compartment  fire  size  or  the  duct  velocity.  In  a 
compartment  without  a  fire  Vf=0.1  m/s  and  with  a  fire  the  velocity  is  calculated  as: 
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where  Azi  is  the  height  of  the  compartment  and  Lj  is  either  one  half  the  diagonal  distance  from 
the  center  of  the  upper  surface  to  its  comer,  or  one  half  a  user-specified  separation  distance. 

4.2.5  UV/IR  Detector 

The  current  model  for  a  UV  or  IR  detector  serves  as  a  placeholder.  In  its  current 
implementation  the  detector  will  alarm  when  the  fire  reaches  a  user  defined  size. 

4.3  Fire  Suppression  Models 
4.3.1  Mist 

Mist  systems  are  assumed  to  extinguish  a  fire  via  displacement  of  oxygen  due  to 
evaporation.  The  mist  model  tracks  both  a  mist  droplet  quantity  and  an  evaporation  quantity  for 
the  compartment  in  which  the  mist  system  is  located.  Currently,  evaporation  of  mist  from 
remote  compartments,  convective  heating  of  droplets,  droplet-Aplet  interactions,  and 
droplet-surface  heat  transfer  are  not  evaluated.  The  mist  model  involves  two  primary 
computations.  The  first  is  the  time  rate  of  change  of  droplet  concentration  resulting  from  mist 
nozzle  discharge  and  droplet  fallout.  The  second  is  the  rate  of  evaporation,  which  is  determined 
by  available  energy  for  evaporation  and  the  compartment’s  equilibrium  vapor  mass  fraction. 

The  droplet  concentration  at  the  end  of  a  time  step  is  given  by: 

-  Mmou.  -  (64) 


Fallout  results  from  mist  droplets  impacting  one  of  the  surfaces  of  the  compartment.  If  the 
droplets  are  traveling  at  their  terminal  velocity,  the  rate  at  which  they  fallout  is  given  by  the  time 
it  takes  the  droplets  to  travel  from  their  location  to  a  surface.  Since  the  mist  spray  and 
compartment  convection  will  distribute  the  mist  widely  throughout  a  compartment,  there  is  no 
fixed  droplet  location;  instead  a  mean  beam  length.  Equation  44,  is  assumed.  The  terminal 
velocity  and  the  mean  beam  length  jdeld  a  time  constant  for  the  droplet  fallout.  Droplet  terminal 
velocity  is  given  by  solving  the  equation  below  for  the  droplet  terminal  velocity,  Va: 

Mdg  =  ^PiCD7rrd^v^,  (65) 


where  Ma  is  the  droplet  mass,  Cd  is  the  droplet  drag  coefficient,  ra  is  the  droplet  radius,  and  p\  is 
the  compartment  gas  density.  The  drag  coefficient  is  a  fimction  of  the  droplet  Reynolds  number 
and  is  given  by  (Clift,  Grace,  and  Weber,  1978): 
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Since  the  droplet  Reynolds  number,  Rea,  is  a  fimction  of  the  droplet  velocity,  the  terminal 
velocity  is  determined  by  iteration.  With  terminal  velocity  known  the  droplet  fallout  rate  is 
given  by: 
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The  evaporation  of  mist  droplets  is  limited  by  four  factors:  the  compartment  equilibrium 
vapor  mass  fraction,  the  energy  available  for  evaporation,  the  mass  of  available  mist,  and  the 
convective/diffusive  limitation  on  the  droplet  vaporization.  The  equilibrium  vapor  mass  fraction 
can  be  given  by  the  Clausius-Clapeyron  equation  (Wylen  and  Sonntag,1986)  where  Xy  is  the 
equilibrium  mole  fraction  for  water  vapor: 
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where  hy  is  the  heat  of  vaporization  of  water,  mw^^o  the  molecular  weight  of  water,  and  Tb  is 

the  boiling  temperatme  of  water.  Assuming  the  remaining  gas  is  air,  this  is  equivalent  to  a  mass 
fraction  of: 
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The  energy  available  to  evaporate  the  droplet  is  the  sum  of  the  droplet  enthalpy  and  the 
convective  heat  transfer  to  the  ^op.  The  convective  heat  transfer  is  determined  from  the 
following  correlation  for  heat  transfer  coefficient  (Cheremisinoff,  1986): 
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The  rate  at  which  the  water  can  evaporate  is  limited  by  the  diffusion  and  convection  of  water 
vapor  away  from  the  drops  surface  as  it  evaporates  (Cheremisinoff,  1986): 
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where  Djj^q  is  the  diffusion  coefficient  of  water  vapor  into  air  (0.256x10’'*  m^/s),  n<i  is  the 

number  of  droplets  available  for  evaporation,  and  Sh<]  is  the  droplet  Sherwood  number  which  is 
given  by: 


Shj  —  2  +  0.6RejSCgj,. , 


where  Sc  is  the  Schmidt  munber. 


(72) 


The  actual  mist  evaporation  is  calculated  iteratively.  Equations  70  and  71  are  used  along 
with  the  total  available  water  mist,  after  accounting  for  any  settling  with  Equation  67,  to 
determine  the  maximum  possible  mist  evaporation.  A  guess  at  the  new  end  of  time  step 
compartment  mass  and  energy  content  is  made  using  the  previous  time  step  mass  and  energy  and 
prior  iteration  change  in  mass  and  energy  along  with  any  change  due  to  the  predicted 
evaporation.  From  this  a  new  compartment  temperature  can  be  determined.  The  equilibrium 
vapor  mass  fraction  for  this  new  temperature  is  compared  with  the  guess  of  the  new 
compartment  vapor  mass  fraction.  If  the  equilibrium  is  exceeded  a  bisection  method  is  used  to 
guess  a  mass  of  evaporation.  This  process  is  repeated  at  each  inner  iteration.  In  this  manner  the 
calculated  evaporation  should  avoid  any  large  violations  of  the  equilibrium  condition. 

4.3.2  Sprinkler 

The  performance  of  a  standard  sprinkler  system  is  evaluated  using  the  method  of  Evans 
(1993).  This  method,  resulting  from  analysis  of  experimental  data,  assumes  that  a  fire  is 
suppressed  by  a  standard  sprinkler  system  with  a  time  constant,  t,  given  by: 


T  =  3w 


w-1.85 


(73) 
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where  w"  is  the  spray  density  in  kg/m^-s.  Suppression  is  calculated  as: 


Qf, suppressed 


(74) 


where  tact  is  the  activation  time  of  the  sprinkler  system. 

4.3.3  Gaseous 

A  gaseous  suppression  system  will  extinguish  a  fire  once  it  displaces  sufficient  oxygen  to 
bring  the  compartment  oxygen  mass  firaction  below  a  specified  lower  limit  for  combustion. 

4.3.4  Hydroflourocarbon  (HFC) 

An  HFC  system  will  extinguish  a  fire  once  the  agent  reaches  a  specified  extinguishment 
mass  fi’action. 

4.3.5  Handline 

No  model  implemented  at  this  time. 

4.3.6  Foam 

The  foam  suppression  model  calculates  a  suppression  factor  which  is  equal  to  the  firaction  of 
the  floor  area  that  is  covered  by  the  foam.  The  area  covered  is  calculated  by  using  the  foam  flow 
rate,  foam  density,  and  foam  layer  thickness  combined  with  the  net  time  the  foam  system  has 
been  discharging.  This  model  serves  as  a  placeholder  for  future  development  of  a  more 
physically  based  model. 

4.3.7  Dry  Chemical 

A  dry  chemical  system  will  extinguish  a  fire  once  the  agent  reaches  a  specified 
extinguishment  mass  fraction. 

4.3.8  Boundary  Cooling 

A  functioning  boundary  cooling  system  is  assumed  to  coat  the  entire  surface  being  cooled 
with  a  film  of  water  of  constant  temperature.  The  boundary  cooling  operates  by  replacing  the 
heat  transfer  coefficient  in  Equations  30  through  32  with  a  constant  heat  transfer  coefficient. 

This  coefficient  is  either  3000  W/m^  or  1000  W/m^  with  the  latter  representing  a  Leidenffost,  or 
film  boiling,  condition  which  is  reached  when  the  surface  temperature  is  more  than  200  K  greater 
than  the  boiling  point. 
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4.4  Fan  Model 


FSSIM  has  three  models  for  evaluating  the  performance  of  a  ventilation  fan^lower.  These 
are  a  constant  pressure  model,  a  quadratic  model,  and  a  user  specified  flow  vs.  pressure  drop 
(head  loss)  model.  For  the  constant  pressure  model,  the  user  specifies  a  fixed  pressure  rise 
across  the  fan.  This  pressure  is  added  to  the  momentum  equations  as  a  source  term,  the  APj  term 
in  equation  4.  For  the  quadratic  model  the  user  specifies  the  maximum  volumetric  flow,  'Vf^^ ,  at 
zero  pressure  drop  across  the  fan  and  the  pressure  drop  at  which  the  flow  rate  becomes  zero, 

^^shutoff  • 


(75) 

At  the  beginning  of  each  time  step  the  volumetric  flow  through  the  fan  is  calculated.  The 
fan  pressure  is  solved  for  and  used  as  the  source  term.  To  prevent  oscillations  in  the  fan  source 
term  the  added  pressure  is  relaxed  between  time  steps  using: 

AP;=pAP^.  +  (l-p)AP;-',  (76) 


APj  -  APj^„^jf 


.  V...^  , 


where  p  is  a  relaxation  factor,  p=0.3. 

For  the  user  specified  model  the  user  inputs  an  array  of  pressure  drop  vs.  volumetric  flow 
rate.  At  each  time  step  the  current  volumetric  flow  through  the  fan  is  determined  fi-om  the  prior 
time  step  velocity  in  the  duct.  A  linear  interpolation  is  performed  to  determine  the  next  time  step 
pressure  which  is  then  relaxed  using  Equation  76. 
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NOMENCLATURE 


5.0 


5.1  Roman 

A  Area  (m^) 

Cd  Discharge  coefficient  or  droplet  drag  coefficient 

c  Specific  heat  (J/kg-K) 

Cp  Constant  pressure  specific  heat  (J/kg-K) 

D  Diffusion  coefficient  (m^/s) 

d  Diameter 

E  Energy  (J) 

F  Function  or  view  factor 

f  Friction  factor 

g  Gravitational  acceleration  (9.80665  m/s^) 

h  Enthalpy  (J/kg)  or  heat  transfer  coefficient  (W/m^-K) 

K  Form  loss  factor 

k  Thermal  conductivity  (W/m-K) 

L  Length  (m) 

M  Mass  (kg) 

mw  Molecular  weight  (kg/mol) 

OD  Optical  depth  (m' j 

P  Pressure  (Pa) 

Pr  Prandtl  number 

q  Heat  transfer  rate  (W) 

R  Real  gas  constant  (8.3 14472  N-m/mol-K) 

Re  Reynolds  number 

r  Radius  (m) 

RTI  Response  Time  Index  (m'^s*^^) 

Sc  Schmidt  number 

Sh  Sherwood  number 

T  Temperature  (K) 

t  Time  (s) 

V  Volume  (m^) 

V  Velocity  (m/s) 

w  Suppression  mass  (kg). 

X  Mole  fraction  (mol/mol)  or  x/z 

X  Position  (m) 

V  Mass  fraction  (kg/kg)  or  y/z 

z  Vertical  position  (m) 

5.2  Greek 

a  Absorptivity  (m'*) 

P  Miscellaneous  multiplier  (i.e.  relaxation  factor) 

A  Change  in 

e  Emissivity  or  wall  roughness  (m) 


34 


p  Density  (kg/m^) 

CT  Stefan-Boltzman  constant  (5 .6704x1  O'*  W/m^-K'*)  or  direction  function  (1  or  -1) 

T  Transmission  factor  or  time  constant  (s) 

Xf  Radiative  fraction 

5.3  Superscripts 

Per  unit  length  (m"') 

Per  unit  area  (m'^) 

Per  unit  volume  (m'*) 

in  Quantity  is  flowing  inward  to  a  node  or  compartment 

n  Next  time  step 

n+  Next  time  step  guess 

n-  Prior  iteration 

n-1  Prior  time  step 

out  Quantity  is  flowing  outward  from  node  or  compartment 


5.4  Subscripts 


act 

air 

b 

c 

d 

eq 

equil 

evaporation 

ex 

f 

fallout 

flood 

full 

gas 

H2O 

i 

j 

k 

1 

loss 

mist 

n 

net 

nozzle 

outflow 

pyro 


Activation 

Air 

Boiling 

Convection 

Droplet,  duct, or  diameter 
Equivalent 

Evaluated  at  liquid- vapor  equilibrium 
Results  from  evaporation  of  a  liquid 
Exchange 
Fire 

Fallout  of  a  liquid  or  solid  aerosol  from  the  gas  phase 

Evaluated  at  the  flooding  criteria  for  vertical  flow  through  an  opening 

Designates  flowrate  for  fan  or  blower  at  zero  head 

Gas 

Water 

Compartment 

Junction 

Compartment 

Liquid 

Mist  droplet 
Duct  node 

Sum  of  incoming  and  outgoing  parameter 
Results  from  injection  through  a  nozzle 
Results  from  flow  out  of  a  compartment  or  node 
Pyrolysis 
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r 


r 

s 

shutoff 

supp 

suppressed 

V 

w 


Radiation 

Soot 

Designates  head  at  which  fan/blower  flow  rate  goes  to  zero 

Suppression  system  parameter 

Suppressed 

Vapor 

Surface 


5.5  Overscripts 

Time  derivative  (s’*) 

~  Linearized  value  at  next  time  step  or  modified  net  radiation  term 
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